Study of two-dimensional electron systems in the renormalized-ring-diagram 

approximation 
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With a super-high-efficient numerical algorithm, we are able to self-consistently calculate the 
Green's function in the renormalized-ring-diagram approximation for a two-dimensional electron 
system with long-range Coulomb interactions. The obtained ground-state energy is found to be in 
excellent agreement with that of the Monte Carlo simulation. The numerical results of the self- 
energy, the effective mass, the distribution function, and the renormalization factor of the Green's 
function for the coupling constants in the range < r 3 < 30 are also presented. 
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Two-dimensional electron systems (2DES) with long- 
range Coulomb interactions, realized in the semiconduc- 
tor heterostructures and inversion layers Q] , continuously 
attract a lot of attentions despite it has been studied for 
more than three decades. Most of the theoretical calcula- 
tions so far were based on the random-phase approxima- 
tion (rpa) HUSH, 

which is expected to be accurate 
in the weak coupling or high density regime. For taking 
into account the strong-coupling effect beyond the PRA, 
the usual way is to adopt the local-field correction (LFC) 
to the RPA 0, 0, S H EH G3 In all these works, the 
Green's functions appeared in the diagrams of the stan- 
dard RPA as well as the one with the LFC are not renor- 
malized. A crucial condition for a better approximation 
is that the renormalized Green's function should satisfy 
some microscopic conservation laws, otherwise the ap- 
proximation may lead to unphysical consequences |12| . 
Among various approximations, the renormalized-ring- 
diagram approximation (RRDA), in which the Green's 
function needs to be self-consistently determined from 
the relevant i nteg ral equations, is well known to satisfy 
this condition |12( and thus should be a sound approach. 
Because tremendous numerical efforts were needed to 
solve the integral equations in RRDA, the physics in this 
approach had only been studied for small coupling con- 
stants EH i an d there existed no solutions to RRDA from 
the intermediate to strong-coupling regimes. To develop 
a numerical scheme to solve RRDA equations for the pur- 
pose of understanding 2DES in a wider range of couplings 
is still a challenging problem in modern many-particle 
physics. In this paper, we use a super-high-efficient algo- 
rithm to solve the RRDA equations. With the solutions, 
the ground-state energy is obtained, and it is in excel- 
lent agreement with the result of the fixed-node-diffusion 
Monte Carlo (MC) simulation E3- ^ n addition, the self- 
energy of the Green's function, the effective mass, the 
single particle distribution function as well as the renor- 
malization factor of the Green's function are also calcu- 
lated. 

We consider a 2DES with density of n at temperature 



T embedded in an uniform neutralizing background of 
positive charge. Through out this paper, we will use 
the units in which Ti — ks = m = a = 1 (with m the 
mass of the electron, and a the Wigner-Seitz radius of 
an electron in the 2DES). The system is characterized by 
two dimensionless parameters: — T/Ep and r s — a/as 
where Ep = nn (= 1 in our units) is the Fermi energy and 
as the Bohr's radius. The electronic Green's function G 
is related to the self-energy E via 



G(k, iw n ) = [iu! n - £r - E(fc, icon)]' 



(1) 



where = — with = k 2 /2 and [i the chemical 
potential, and u„ is the fermionic Matsubara frequency. 
For brevity, we hereafter will use k = (k,iu> n ) for the 
argument unless stated otherwise. By the RRDA [l2j| . 
the self-energy E is given by 



V/3 



q> l-V(q)x(qY 



(2) 



where (3 — l/T, V(q) — 2%r s /q is the Coulomb interac- 
tion, V is the area of the system, and 



Xil) = ypY,G{k)G{k + q) 



(3) 



is the electron polarizability, q = (q,iQ, m ) with Q m the 
bosonic Matsubara frequency. The chemical potential /i 
is determined by the electron density, 



(4) 



By solving Eqs. ^ ne f unc ti° ns G and E, and the 

chemical potential /i can be self-consistently obtained. 

In the numerical procedure, the most time-consuming 
computations are calculations of E and x because of the 
summation over the Matsubara frequency and the in- 
tegration over the momentum. We have developed an 
effective algorithm for doing such a computation. In 
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FIG. 1: (Color online) Self-energy E(fe, iu) n ) as function of n 
a,t k — kF, — 0.03, and r a = 10. 



this method, the summation is taken only over some 
selected Matsubara frequencies distributed in L succes- 
sively connected blocks each of them containing M fre- 
quencies, with each term under the summation multiplied 
by a weighting factor. For details of this algorithm, the 
reader is referred to Ref. In the present calcula- 

tion, we have used the parameters [h,L,M] = [2,15,5] 
for selection of the Matsubara frequencies, where h is 
the integer-parameter that the stride in the £-th block 
is h^ -1 '. The total number of the selected frequencies 
is L(M — 1) + 1 = 61 with the largest number N ~ 
2 L (M - 1) = 2 17 for the cutoff frequencies fl N = 2NttT 
and ojn = (2N — l)nT. We will see later that the 61 
Matsubara frequencies are sufficient to describe the self- 
energy and thereby the Green's function. Instead of sum- 
ming over N Matsubara frequencies by a usual method, 
we here take the summation only over the 61 ones. There- 
fore, the efficiency of the algorithm for the present calcu- 
lation is N/61 ~ 2 11 ! For the lowest temperature consid- 
ered here, 9 — 0.03, we have uj^/Ef ~ 2.47 x 10 4 . In our 
calculation, the largest is cm = k 2 M /2 with k,M — 50. 
Therefore, the cutoff ujn/^m ~ 20 is sufficiently large. 

On the other hand, the momentum-space convolution 
integrals in Eqs. J5J) and can be easily performed by 
Fourier transforming into real space. For illustrating the 
two-dimensional Fourier transform, we here discuss the 
transformation of the Green's function. In the real space, 
the Green's function is given by 



G(r, iu> r , 



dk 
2^ 



kG(k, iuj n )Jo(kr) 



(5) 



where Jo is the first kind Bessel function of order zero. 
However, Eq. J^Jl is not in a favorite form for the numer- 
ical computation since kG(k 1 iu> n ) behaves like 0(1/ k) 
at k — > oo and jo(kr) is an oscillatory function; the 
integrand decreases so slowly at large k that a precise 
numerical result is hardly obtained. By observing the 
asymptotic behavior, 



with G°(k, iu> n ) — l/(iu) n +/j,Q — k 2 /2) and respectively 
the Green's function and the chemical potential of the 
free electrons, we choose G° as the auxiliary function for 
G. The Fourier transform of G (k,iui n ) is given by 



G°(r,iiu n ) = -K (pr)/ir 



(7) 



where Kq is the second kind modified Bessel function of 
order zero, and p = (— 2/io — i2oj n ) 1 / 2 . Thus, equation 
(JSJ can be rewritten as 

f°° dk 

G{r,iu) n )=J —kSG(k,iuj n )J (kr)-K (pr)/TT (8) 

with 8G(k, iu> n ) — G(k, iw n ) — G°(k, iu) n ). Now, the inte- 
grand in Eq. JHJ drops fastly with kSG(k, iuj n ) ~ 0(k~ 3 ) 
at k — > oo. Equation (JSJ can be further reformed for deal- 
ing with the fast oscillatory behavior of Jo(kr) at large 
r. For doing this, we consider the integral 



(9) 



g(r) = / dkf(k)Mkr). 
Jo 

For Jo(fcr), we choose the auxiliary function as 



A(kr) = 



kr + c 
kr 



8(kr + c) 



{sin(fcr) + cos(fcr) 

sin(fcr) - cos(fcr)]} (10) 



with c = 1/it. A(kr) has the same asymptotic behavior 
as Jo(fcr) at kr — > oo. In addition, A(0) — Jo(0) = 1. 
Equation can be rewritten as 



9(r) 



dkf(k)[J (kr)-A(kr)} 



f(k) 



f(k) 



c I dk 

o yfcr + c 



kr 



i c I dk , i , 



ki 



8(kr + c) 



cos(fcr) 
2 j sin(fcr). 



G(k,iuj n ) — > G°(k, iu n ), 



at k 



oo or n 



(6) 



Now, the first integral in the above equation can be in- 
tegrated by the simple numerical method since the fast 
oscillatory behavior of Jo is considerably cancelled by 
A(kr). On the other hand, the second and third Fourier 
integrals can be numerically integrated using the Filon's 
method. In our calculation, the momentum integration 
is over the range < k < 50. Because the electron 
distribution function varies drastically around the Fermi 
momentum kp = \/2, we divide the fc-range into three 
segments: [0,/cf — &o], [ fcjr — fcp, fcjr + fcp], and [fcjr + fcp, 50] 
with ko = v°0/(l + a/ 326*) . The fc-integrals are calcu- 
lated with 200 meshes in each segment. In the real space, 
r is confined to < r < 80, which is divided into three 
segments: [0,5], [5,30], and [30,80], with 300 meshes in 
each segment. 

In the above example, the accuracy and efficiency of 
numerical integration is improved by choosing the proper 
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FIG. 2: (Color online) Ground-state energy in atomic unit 
(a.u.) as function of r s . The present result (RRDA) is com- 
pared with that of the RPA, GB 0, and MC [3]. 



FIG. 3: (Color online) Effective mass m*/m as function of r 3 . 
The present result at 9 = 0.03 is compared with that at 9 — 
obtained in Ref. [3 (FS) with a different numerical method. 



auxiliary function. Analogously, we can choose the auxil- 
iary function for the summation over the Matsubara fre- 
quency as well. For example, the auxiliary function for 
calculating x(l) is chosen as G° (k)G° (k + q) that is the 
asymptotic limit of G(k)G(k + q) at w„ — * oo or k — > oo. 
On the other hand, the free-particle polarizability x (<?) 
can be easily obtained by 



X°(q,i£l m ) = -r- 

4:TT 



de- 



Y{q, n r , 



1 



cosh%- Mo )/2T] 



(11) 



with Y(q,Q m ) = (^x 2 + y 2 + x) 1 / 2 and 1 = 1/2- 
n 2 n /2e 2 q - 2e/e q , y = a m /e q . 

By iteration, we have solved Eqs. Q-Q for < r s < 
30 at finite temperatures. Firstly, in order to illustrate 
the validity of our selection for the Matsubara frequen- 
cies, we show in Fig. ^ the result for the self-energy 
E(fc, iw n ) for k = hp, 9 — 0.03, and r s — 10. The sym- 
bols represent the values of E(fc, iw n ) calculated at the 
selected n's. This is a typical result of T,(k,ioj n ) in the 
strong-coupling regime. As seen from Fig. ^ the se- 
lected frequencies are sufficient for describing E(&,mj„); 
though the selected uj n 's are sparsely distributed at large 
n, meanwhile, the function E(fc, iu> n ) varies slowly. At 
uj n — > oo, the imaginary part of E(fc,iu; n ) vanishes, but 
the real part remains finite. Actually, at oj n — > oo, we 
get Yi(k,iu) n ) — > the exchange part of the self-energy. 

In Fig. [3 the ground-state energy per electron in 
atomic unit is shown as a function of r s . The results 
of present calculation are obtained by extrapolation of 
the finite-temperature values. Our main result denoted 
by RRDA is compared with the RPA, the high-density 
expansion |3J by the scheme of Gcll-Mann and Brucckner 
(GB) [HI, and the MC simulation Q. In the weak cou- 
pling regime, both the RRDA and RPA reproduce the 
GB result very well. However, for strong coupling, the 
RRDA is much closer to the MC than the RPA. 

Shown in Fig. [21 is the effective mass to* /to as a func- 
tion of r s at 9 = 0.03. The values calculated by Faleev 



and Stockman (FS) with a different numerical scheme 
for T = is also depicted for comparison. The present 
result shows that to* /m is a monotonically decreasing 
function of r s in a wide range of r s . This is different 
from the RPA (not shown here) by which to* /to de- 
creases at small r s < 0.1, but then increases with r s 0. 
The present result on to* /to seems to be contradictory 
to the explanation to the experiments on the tempera- 
ture dependence of the Shubnikov-de Haas oscillations in 
a 2DES |ri|. The experiments were explained as to* /to 
increasing with r s . To this point, there may be two pos- 
sibilities: (1) the RRDA may not be good enough, or 
(2) the experimental results may require a different in- 
terpretation; for example, the conductivity formula to fit 
the experiments was based on non-interacting electrons 
except the mass is assumed to be renormalized. Further 
investigation is needed to get a clear understanding of 
this problem. 

The effective mass to* /to has been extensively inves- 
tigated by many calculations using the RPA. With the 
RPA self-energy, to* /to is calculated by the on-shell 
scheme 0), which is not a self-consistent approxima- 
tion for the quasiparticle energy. Contrary to the RPA, 
the RRDA is a conservin g ap proximation for the single- 
particle Green's function |l21: the Green's function sat- 
isfies the Luttinger theorem |l8l fl(| by which the Fermi 
surface is unchanged for the 2DES. To see this, we show 
in Fig. 0] the retarded self-energy ~E r (k,E) as function 
of E at k = k F for 6 = 0.03 and r s = 10. E r (k,E) 
is obtained by the Pade approximation j^. From Fig. 
Eta), it is seen that ReE r (fcF, E) is essentially a lin- 
ear function of E in the neighborhood of E = 0. It 
is clear that ReE r (fcf,0) + £fc F = 0, which means that 
the quasiparticle energy E vanishes at the Fermi sur- 
face consistent with the Luttinger theorem. The equal- 
ity ReE r (fc^,0) + £fc F = originates from the fact that 
ReE r (fci?,0) is compensated by the shift of the renor- 
malized fi from /xq. The imaginary part shown in Fig. 
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FIG. 4: (Color online) Self-energy Y, r (k,E) as function of £ 
at k = k F for 6 = 0.03 and r s = 10. (a) ReT, r (k Fl E) + £ kF 
and (b) the imaginary part. 

Efb) has a parabola form with a small value at E = 0. 
At T = 0, ImS r (fcF,-B) cx -E 12 at small E 18]. The 
small negative value of lm£ r (fci?,0) shown in Fig. Efb) 
is a finite-temperature effect. (All these facts indicate 
that our numerical computation is very accurately per- 
formed.) On the other hand, since the bare chemical 
potential /xq is used in the RPA calculation, ReT, r (kp, 0) 
cannot be compensated and thereby the quasiparticle en- 
ergy does not vanishes at hp. 

In Fig. |5j we exhibit the result for the distribution 
function n{k) at = 0.03 for various r s . It is well known 
that n(k) at r s — is the Fermi-Dirac distribution func- 
tion F (k) for the free electrons. At very small r S) n(k) 
is close to F°(k). With increasing r s , n(k) becomes very 
different from F (k); it is gradually suppressed at k < kp 
while growing up at k > kp. At the Fermi momentum, 
n(k) decreases dramatically, but seems continuously. At 
zero temperature, n(k) should have a discontinuity at 
k = kp. The continuous behavior appeared in Fig. [S] 
is attributed to the finite-temperature effect by which 
the discontinuity is rounded. Shown in the inset is the 
renormalization factor Z of the Green's function. Z is 
associated with the abrupt drop in n(k) of the ground 
state at kp. Within the range of r s studied here, the 
magnitude of Z is consistent with the sharp drop in n(k) 
at kp. For comparison, the result by FS is also shown 
in the inset. By their method they could solve the 
RRDA equations only within r s < 2.62. Clearly, their 
result is reproduced by the present calculation. 

In summary, we have solved the RRDA equations for 
2DES at r s < 30 with our super-high-efficient numerical 
algorithm. The obtained ground-state energy is in ex- 
cellent agreement with the Monte Carlo result. For the 
effective mass m*/m, the RRDA result shows m* /m is 
a decreasing function of r s , which is considerably differ- 
ent from the RPA. The distribution function n{k) and 
the renormalization factor Z of the Green's function are 
calculated for a wide range of r s . 
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sion on the related problems. This work was supported 
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FIG. 5: (Color online) Distribution function n(k) at 9 — 0.03 
and various r s . The inset is the Z- factor as function of r s . 
The circles represent the results of the present calculation, 
and the squares are that of Ref. Fill 
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